#!/bin/bash
set -e

while getopts  "s:p:t:" opts
do
        case  $opts  in
        s)
			sample_name=$OPTARG
			;;
		p)
			out_prefix=$OPTARG
			;;
		t)
			total_reads=$OPTARG
			;;
		\?)
			echo `basename $0` '[-p out_prefix] [-t total_reads] <in.bam>'
			exit 1
			;;
		
        esac
done
shift $(($OPTIND - 1))


if [ $# -lt 1 ]; then
        echo `basename $0` '[-p out_prefix] [-t total_reads] <in.bam>'
        exit 1
fi

# if [ -z $sample_name ]; then
        # sample_name=sap
        
# fi

if [ -z $out_prefix ]; then
	out_prefix=1
fi

#-----------------------------------------------
#-----------------------------------------------
. /mnt/ilustre/app/medical/tools/.var #---------
#-----------------------------------------------
#-----------------------------------------------

log=.log
if [ ! -e "$log" ]; then
	:> $log
fi

reads_len=135
genome_len=3095677412



# trc total reads counts
# samtools view -c $1 > $out_prefix.trc.txt

if test -z $total_reads
then
total_reads=$(samtools view -c $1)
fi
echo &>>$log
echo &>>$log
echo total reads of $1: $total_reads &>>$log

((depth=${total_reads}*reads_len/genome_len))

echo &>>$log
echo &>>$log
echo $1 depth is: $depth &>>$log

if test $depth -gt 1
then

#---------------------------------------

bedfile=b37.chosen_chr.bed
bedfile1=b37.rand.site.chosen_chr.bed
chosen_chr=(7 13 X)
for chr in ${chosen_chr[@]}
do
cat ${data_path}/intervals/1/b37.chr1y.bed |awk '$1~/^'"$chr"'$/' >> $bedfile
cat ${data_path}/intervals/1/cnv_sc/b37.rand.site.bed |awk '$1~/^'"$chr"'$/' >> $bedfile1
done

chosen_chr_len=$(cat $bedfile |awk 'BEGIN{sum=0}{sum+=$3}END{print sum}')

#---------------------------------------


for i in $(seq $[depth-1])
do

echo &>>$log
echo &>>$log
echo ${i}X &>>$log
# ((sample_fraction=$i/$depth))

sample_fraction=`echo "scale=8;$i/$depth"|bc` # cov_dep_sub_sam.sh is wrong in this

# bamtools random -in $1 -out $out_prefix.${i}X.bam -n $extract_num 2>>$log
samtools view -@ 8 -b -s $sample_fraction -L $bedfile $1 > $out_prefix.${i}x.sub.bam 2>>$log
samtools index $out_prefix.${i}x.sub.bam 2>>$log

# samtools sort -@ 9 -m 3G -T $out_prefix -O bam $out_prefix.${i}X.bam \
# > $out_prefix.${i}X.sort.bam \
# 2>>$log

# samtools index $out_prefix.${i}X.sort.bam 2>>$log

# rm $out_prefix.${i}X.bam

summ_sc_gatk.sh -p $out_prefix.${i}x.sub $out_prefix.${i}x.sub.bam $bedfile1



done

mv $bedfile1 b37.rand.site.chosen_chrs.bed
mv $bedfile b37.chosen_chrs.bed

else
echo &>>$log
echo &>>$log
echo Initial depth is not more than 1X, do nothing. &>>$log
fi

echo &>>$log
echo &>>$log
echo done &>>$log

exit 0








